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Binary encounters between spherical particles in shear flow are studied for a system 
bounded by a single planar wall or two parallel planar walls under creeping flow con- 
ditions. We show that wall proximity gives rise to a new class of binary trajectories 
resulting in cross-streamline migration of the particles. The spheres on these new trajec- 
tories do not pass each other (as they would in free space) but instead they swap their 
cross-streamline positions. To determine the significance of the wall-induced particle mi- 
gration, we have evaluated the hydrodynamic self-diffusion coefficient associated with a 
sequence of uncorrelated particle displacements due to binary particle encounters. The 
results of our calculations quantitatively agree with the experimental value obtained by 
Zarraga fc Leig hton (2002) for the sclf-diffusivity in a dilute suspension of spheres under- 



going shear flow in a Couette device. We thus show that the wall-induced cross-streamline 
particle migration is the source of the anomalously large self-diffusivity revealed by their 
experiments. 



1. Introduction 

Random displacements resulting from particle encounters in suspension flows lead to 
hydrodynamically induced particle migration with respect to the local suspension velocity 
(jEckstein et aLlll977tlLeighton fc Acrivosll98~7llBossis fc Bradvll98"7h . For non-Brownian 



particles under creeping- flow conditions, hydrodynamically induced migration constitutes 
an important mechanism of particle redistribution in the suspending fluid. Thus studies 
of particle encounters in flowing suspensions are vital both from the fundamental and 

practical point s of view. 

As shown bv lBatchelor fc Green two spheres passing each other in unbounded 



shear flow return to the initial transverse (cross-streamline) positions after a binary en- 
counter is completed. This behavior follows from the flow reflection symmetry of Stokes 
equations and a reflection symmetry of the system. Since on open binary trajectories 
there are no cross-streamline particle displacements, for perfect spheres in free space the 
transver se hydrodynamic-diffusion process requires encounters of at least three particles 
see, eg. lAcrivos. Batchelor. Hinch. Koch fc Mauril f 992 ; Wang. Mauri fc Acrivod 19961 . 



1998; iDrazer. Koplik. Khusid fc Acrivosll2002f h It follows that the self-diffusion coefficient 
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scales as D s ~ 0(</> 2 ) in the low-concentration regime, and the 0(<p) contribution asso- 
ciated with binary encounters arises only in the presence of non-hydrodynamic forces 
that remove the flow-reflection symmetry of the problem. For non-Brownian particles 
non-hydrodynamic interactions often result from direct particle contacts due to surface 
roughness. 

Several years ago IZarraga fc Leightonl (|2002l ) measured the 0{<f) contribution to the 
shear-induced self-diffusivity D s for a dilute suspension of spheres undergoing shear flow 

in a Couette device. The experiments yielded a surprising result: the self-d iffusion coeffi- 

cient was nearly an order of mag nitude higher than the theoretical estimate (|da Cunha fe Hinch 
1996; IZarragafc Leightonl [20011) for rough spheres with the roughness amplitude corre- 



sponding to the experimental system. Several possible causes of the anomalous self- 
diffusivity (such as inertial lift and non-Newtonian effects) were examined, but none of 
them was sufficient to explain the anomaly. 

The suspension in the Couette device used in experiments of IZarraga fc Leighton 
( 2002J) was bounded by nearly flat parallel walls separated by a relatively small dis- 
tance H = 20c? (where d is the particle diameter). In the analysis of experimental results 
it was assumed that the effect of walls on the self-diffusivity coefficient must be negligi- 
ble. Since the fore-aft symmetry of the system would ensure that after passing each other 
the particles would return to their original streamlines (as they do in infinite space), it 
seemed unlikely that the walls were the cause of the enhanced self-diffusivity. 

In the present paper we show that hydrodynamic interactions of particle pairs with 
the confining walls result in cross-streamline particle displacements in binary encounters 
in shear flow. Specifically, we have found a new class of binary trajectories: the particles 
on such trajectories initially approach each other, but then they move across the channel 
in opposite directions and separate without passing each other (unlike the particles in 
free space). 

The new class of trajectories results from a sign change of the transverse component 
of the relative particle velocity. Such sign changes in wall-bounded systems we re first 
pointed out in our recent study (jBhattacharva. Blawzdziewicz fc Wainrvbll2005aj ) . Thus 
seemingly subtle features of particle mobilities produce a significant qualitative effect: 
our explicit calculations show that after the new class of particle trajectories is included 
in an estimate for the shear-ind uced self-diffusivity, a quan titative agreement is obtained 
between the measurements of Zarraga fc Leightonl ( 20021 ) and the oretical predictions. 
Hence , the paradox of the unusually large self-diffusivity observed by lZarraga fc Leighton 
( 20021 ) has been resolved. 

The system considered in our paper is defined in section [51 The new class of binary 
trajectories resulting in cross-streamline particle displacements and the physical mech- 
anism leading to this behavior is discussed in section [3l In section |4] we analyze the 
consequences of the new pair trajectories for cross-streamline particle migration in sus- 
pensions of spheres in the low-concentration regime. Our findings are summarized in 
section [SJ 



2. Definition of the system 

We consider the dynamics of binary encounters of spherical particles undergoing sta- 
tionary shear flow in a space bounded by a single planar wall or two parallel planar 
walls separated by the distance H . We focus on configurations where the particle-wall 
separation is comparable to the sphere diameter d. The suspending fluid has viscosity 77, 
and creeping flow conditions are assumed. 

We use a coordinate system where the walls are in the x-y planes. The positions of 
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Figure 1. System definition and schematic representation of a non-swapping binary 

trajectory. 



particle centers are denoted by and the linear and angular velocities of the particles 
by Uj and fij, respectively, where i = 1,2. The walls are at z — and z = H . Particle 
encounters are described using the relative coordinates Ar = (Ax, Ay, Az), where 



Ar = Y2 — ri 

is the relative position vector centered on particle 1. 
The unperturbed fluid velocity 

v cxt (r) = jze x 



(2.1) 



(2.2) 



(where 7 denotes the shear rate) points in the lateral direction x, and varies in the normal 
direction z. The flow occurs due to the motion of the upper wall with velocity 



jHe x 



(2.3) 



Particles are torque and force free. No-slip boundary conditions are imposed at both 
walls and at the particle surfaces. 

In our calculations of particle trajectories, interparticle and particle-wall hydrody- 
namic interactions are accurately taken into account. For a one- wall system the particle 
velocities are evaluated using a reflection technique ( Cichocki et al. 2000), and for the 



two- wa ll system we use the Cartes ia n- repre sentation method, recently developed by our 
group ( Bhattacharva et al. 2005 6l fa. 2006a ). The equations of motion = (where 
the dot denotes the ti me derivative) are integrated using a Runge-Kutta algorithm with 
an adaptive time step (jPress et aZ.lll992n . 

The geometry of the system and a sketch of trajectories for two spheres that pass 
each other are shown in figure [TJ The trajectories in this figure are out of scale, and an 
accurate representation of the trajectories can be seen in figure El 



3. Effect of walls on the dynamics of binary collisions 

3.1. The morphology of open trajectories 

An analysis (|Lin. Lee fc S at herl [197(1 iBatchelor fc Greed fl972l IZinchenko|[l984l ) of the 
relative motion of a pair of spheres in unbounded shear flow (|2.2p indicates that all open 
trajectories start from Ax = —00 and extend to Ax — +00 (in a reference frame fixed on 
the initially slower particle 1). During the approach of the particles their relative vertical 
offset Az increases, reaching a maximum when the particle pair crosses the vorticity- 
gradient plane Aa; = 0. However, after the particles have separated, they return to their 
initial cross-streamline positions (y, z), in accordance with the flow-reversal symmetry of 
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Figure 2. Examples of non-swapping (top panel) and swapping (bottom panel) trajectories 
for two spheres in shear flow (|2.2[) between parallel planar walls. The walls are at z — and 
z = 5d, where d is the particle diameter. The lower wall is at rest, and the spheres are moving 
in the flow-gradient plane x-z. Equal-time positions of sphere centers are connected by dotted 
lines. The inset shows the blowup of the trajectory-intersection region. 



Stokes equations and the symmetry of the system with respect to the reflection of the x 
coordinate. 

For sufficiently large initial vertical offsets Az particle trajectories in wall-bounded 
systems are qualitatively similar to those in free space. The only distinctive feature of 
the trajectory depicted in the top panel of figure is that Az decreases before reaching 
the maximum for the particle pair in the symmetry plane Aa; = 0, while in free space 
Az would increase monotonically. 

However, for smaller initial values of Az we find an entirely different behavior: for the 
new kind of open trajectories shown in the bottom panel of figure[2]the offset Az changes 
sign before the relative position Aa; = has been reached. Accordingly, the component 
AU X of the relative velocity 

AU = U 2 -Ui (3.1) 

also changes sign. The particles do not pass each other - they turn around and separate, 
maintaining Ax < for the whole trajectory. As a result, the spheres do not return to 
their initial streamlines at long times but instead they swap their vertical coordinates 
z. Such particle encounters result, therefore, in displacements of the suspended particles 
across streamlines of the external flow. 

The difference in topology of pair trajectories in unbounded space and in a parallel- 
wall channel is clearly seen in figure [3j where the relative particle motion is depicted 
in the reference frame centered on one of the particles. In free space (top panel) the 
contact surface |Ar| = d is surrounded by an envelope of closed orbits, and all open 
trajectories correspond to particles passing each other. In contrast, in the wall-bounded 
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Figure 3. Relative trajectories for pairs of spheres undergoing evolution in shear flow (12. 2[) in 
free space (top panel) and in a channel of width H/d — 5 (middle panel). The spheres move in 
the flow-gradient plane x-z, and the trajectories are shown in the relative coordinates ()2.f I) for 
different initial vertical offsets of particle centers Az. In the wall-bounded system, the sphere 
with the larger initial value of the coordinate z starts at a distance z/d — 2.4 from the lower wall 
(as in figure Bottom panel shows the evolution of the dimensionless interparticle gap (|3.2|l 
for the trajectories displayed in the middle panel. The gap evolution for swapping trajectories 
is represented by the heavy lines. 



system (middle panel) there also exists a region of swapping trajectories, delimited by 
the critical trajectories crossing the plane Az = at the points where AU Z = 0. 
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Figure 4. Vertical component of relative velocity (13. If) for two spheres aligned along the flow 
direction x, versus particle separation Ax. Different lines correspond to different positions z 
of the particle pair with respect to the lower wall. The left panel shows results for a two-wall 
system with wall separation H/d = 5 and the right panel for a one-wall system. Positive values 
of AUz correspond to the swapping-trajectory domain. 



We emphasize that the swapping trajectories do not violate any symmetries of the sys- 
tem. Since the particles do not pass each other, individual trajectories are not symmetric 
with respect to the reflection x — > —x. However, for each trajectory in the halfspace 
Ax < there is a corresponding reflected trajectory in the halfspace Ax > 0. 

The position-swapping trajectories have several important consequences. First, they 
contribute to cross-streamline particle migration in dilute suspensions bounded by planar 
or nearly planar walls. In particular, this migration mechanism expla ins the enhanced 
hydrodynamic self-difTusivity observed by Zarraga fc Leightonl (|2002T ). as discussed in 
section |U 

Another important consequence of the swapping mechanism is that it prevents near- 
contact particle encounters. The results depicted in the middle panel of figure [3] indicate 
that spheres on the swapping trajectories approach each other less closely than spheres 
on the usual trajectories (i.e. when the particles pass each other). A detailed view of the 
evolution of the dimensionlcss gap between the particles, 



e=|Ar|/d-l, 



(3.2) 



is shown in the bottom panel of figure [3l The results (shown for each of the trajectories 
displayed in the middle panel) indicate that the gap on non-swapping trajectories may 
decrease to about e w 10 -5 , whereas on the swapping trajectories e always remains 
of order one. Thus, the swapping mechanism may prevent particle aggregation in the 
presence of short-range attractive forces. 

3.2. Mechanism for particle swapping 

As shown in the middle panel of figure [3l the swapping trajectories occur when the 
vertical component AU Z of the relative particle velocity (|3.1|) is negative in a subdomain 
of particle configuration with Aa; < and positive in the corresponding subdomain of 
the region Ax > 0. Consider, for example, a pair of spheres located in a flow-vorticity 
plane z — const, and assume that Ax > 0, i.e., sphere 1 is to the left of sphere 2. If 
AU Z > in this configuration, sphere 2 moves from the region below sphere 1 (where 
AU X < 0) to the region above sphere 1 (where AU X > 0). Thus, the particles initially 
approach each other but then they separate without passing each other - they swap their 
vertical positions instead. 
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Figure 5. Fluid streamlines around a single sphere in external shear flow (|2.2[l in free space (top 
panel) and in a channel of width H/d = 5 (bottom panel). Trajectories are depicted in coordinate 
system centered on the sphere. The particle in the channel is in the midplane z = H/2. 

The dependence of the relative vertical velocity AU Z on particle separation Ax is 
depicted in figure |4] for a particle pair aligned in the flow direction. The results indicate 
that in wall-bounded systems the relative vertical particle velocity changes sign at a 
critical particle separation A cr ;t that depends on the channel width and the position of 
the particles with respect to the walls of the channel. For Ax > A cr ;t the sign of AU Z is 
consistent with the topology of the swapping trajectories, and for Ax < A cr j t the sign of 
AU Z corresponds to closed orbits. No such sign change, however, occurs in free space. 

Although particle behavior on swapping trajectories seems at first counterintuitive 
(the vector connecting the particle centers rotates with an angular velocity opposite to 
the vorticity of the external flow) it can be qualitatively explained by considering an 
analogy with shear flow in a channel blocked by an obstacle. If the channel is completely 
blocked, the fluid dragged by a wall towards the obstacle is deflected to the other side 
of the channel and returns along the opposite wall. If the channel is partly blocked (by 
one of the spheres in our system) the streamlines that directly approach the obstacle are 
deflected and exhibit the recirculation behavior. Analogously, in a two-sphere system the 
swapping motion of each particle is produced by the flow deflected by the other particle. 

A more detailed, quantitative explanation of the swapping mechanism is obtained by 
analyzing the effect of the walls on the flow pattern vi around a single sphere in shear 
flow in a weakly confined system. To leading order in the particle-particle and particle- 
wall separations, the fluid velocity vi and the velocity of the second particle are identical. 
Therefore, for weak confinements the comparison of the particle and fluid velocities has 
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Figure 6. Explanation of the position-swapping mechanism: the wall reflection of the stresslet 
flow produced by one of the spheres drives the other sphere across the channel, causing reversal 
of the relative particle motion. 



a quantitative meaning (assuming that the interparticle distance remains sufficiently 
large), and for strongly confined systems, pair trajectories are expected to qualitatively 
resemble the streamlines of vi . 

The flow pattern around a single force- and torque-free sphere in an unbounded and 
in a wall-bounded shear flow is illustrated in figure [5J The streamlines are depicted in a 
the coordinate system centered on the particle. The top panel shows the familiar velocity 
field for a particle in free space, and the bottom panel represents the corresponding flow 
for a sphere in the middle of a channel with wall separation H = 5d. As expected, in 
the wall-bounded system there is a region of recirculating streamlines analogous to the 
swapping trajectories shown in figure [31 

The hydrodynamic mechanism that leads to the recirculating streamlines (and hence 
to the particle position-swapping) can be explained by examining the wall reflection 5vl 
of the perturbation flow <5vi produced by a sphere in external shear flow, as schematically 
represented in figure [6] Only the reflection from the lower wall will be discussed because 
the upper wall produces an analogous effect. In our analysis we assume that the wall is 
in the plane z = and the particle is at r = (0, 0, Zi), where z\jd 3> 1. 

The perturbation flow 5vi that results from the scattering of the external flow l]2.2[) 
by the sphere can be expressed by the formula (|Kim fc Karrilal[l99l 



6vi(t) = -hd 3 



2 f 4 8 f 5 



(3.3) 



where r = (x, y, z) (with z = z — z%) is the position of a field point r with respect to the 
particle center, f = |f|, and f = f/f. The flow 5vi is a superposition of the 0(f~ 2 ) radial 
stresslet contribution (the first term on the right-hand side of the above equation) and 
the 0(r~ A ) potential contribution (the second term). 

The streamlines of the perturbation velocity field (|3.3|) are shown in the left panel of 
figure [7] for a system with particle- wall separation Z\/d = 5. In the plane z = (the 
horizontal plane passing through the particle center), the radial stresslet contribution 
vanishes, and the shorter-range potential-flow contribution is vertical, with an orientation 
consistent with the vorticity of the external flow. 

Relation ()3.3|) indicates that at the wall surface (i.e. the plane z = — z±) the pertur- 
bation flow <5vi points downwards for x > and upwards for x < 0. To ensure that 
there is no fluid flux through the wall, the flow reflected from the wall, iJv^, must point 
in the opposite direction. Therefore, assuming that the vertical component of Svl is a 
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Figure 7. Scattered flow (|3.3[) produced by a sphere in shear flow (|2.2[l (left panel) and reflection 
<5v* of this flow from the wall located at z/d — —5 (right panel). The direction of the flow 8vl 
is consistent with the swapping mechanism illustrated in figure [6] 



monotonic function of z (at least up to the particle position), we find that the reflected 
flow Svl has the direction needed to produce the recirculating streamlines (cf. figure [6]). 

It can be shown that 5vl decays as 0(f' 2 ) with the distance f' from the image 
singularity at r = (0,0, —zi) (the image flow involves the reflection of the stresslet). It 
follows that for z — the flow 8v\ dominates the oppositely directed <3(f -4 ) potential- 
flow term of <5vi, except in the region near the particle. This behavior is consistent with 
the pattern of the streamlines depicted in the bottom panel of figure [5] and with the 
topology of pair trajectories represented in the middle panel of figure 03 

The above essential features of the reflected flow can be explicitly seen for a simple case 
of free interface. To the leading order in z\/d, the flow reflected from such an interface 
is given by the expression 

K(r) = 7d 3 ^?', (3.4) 
lb r' 

where z! = z + Z\ and r = (x, y, z) denote the z coordinate and position vector in the 
reference frame centered at the image singularity. Note that for positions close to the 
particle, the magnitude of the flow (|3.4p is 0(d 2 x/f' 3 ), because both d<f' and x <C f'. 

For rigid walls, the crucial features of the reflected flow (i.e., the direction consistent 
with the swapping particle motion, the overall (d/f') 2 decay, and the small additiona l 
factor 0(x/f') for 5v\ z near the plane x = 0) can be derived using iLorentd l|l907t ) 
reflection formula, as discussed in Appendix [A] These features can also be seen from the 
streamlines of <5v* depicted in the right panel of figure [7] 

In the above analysis of the particle-swapping mechanism we have focused on a sin- 
gle reflection 5w\ of the perturbation flow from a wall. However, complete multi-body 
hydrodynamic particle interactions in wall presence are accurately accounted for in our 
quantitative calculations. While our paper is focused on particle motion in the creeping 
flow regime we would also like to note that fluid reversal zones and position-swapping 
pair trajectories similar to those seen in figures [3] and El were observed at finite Reynolds 
numbers for unconfined con figurations ( Mikulencak k, Morris! 20041 Subramanian fc L 
I2006L lKulkarni-Morrisll2007l ). 
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Figure 8. Location of zeros of the relative vertical particle velocity AU Z and sign of AU Z for 
two spheres in the midplane z = H/2 of a channel with wall separation H/d = 5. Shaded region 
(of diameter d) corresponds to overlapping spheres. Solid lines indicate location of zeros that 
delimit the region of swapping trajectories, and dotted lines represent zeros in the periodic-orbit 
domain. 
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Figure 9. Radius A or it of the circle of zeros of AU Z . Top panel shows the normalized radius 
versus channel width H for a particle pair in the midplane z = H/2. Solid line represents 
numerical results and dashed line the asymptotic scaling (|3.6[) . Bottom panels show A cr i t versus 
position 2 of the particle pair for several channel widths H/d (as labeled). The results in the 
left panel are scaled to emphasize asymptotic behavior (|3.6p . and the results in the right panel 
are shown unsealed to emphasize the near-wall behavior. 
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Figure 10. Upstream domains (in transverse coordinates) corresponding to swapping trajecto- 
ries, for different values of the channel width (as labeled) . The top panel depicts these domains in 
coordinates (Ay, Az) for particles moving symmetrically with respect to the midplane z = H/2. 
The bottom panels show the relative vertical offset Az versus the upstream vertical position Z\ 
of the particle closer to the wall, for particles moving in the flow-gradient plane Ay = 0. The 
results in the top and left bottom panels are rescaled to emphasize the asymptotic behavior 
(|3.7I) . and the results in the bottom right panel are shown unsealed to emphasize the near- wall 
behavior. The curves represent the initial coordinates of the limiting trajectories, and the regions 
below the curves correspond to the swapping trajectory domains. 



3.3. The domain of swapping trajectories 

A pair of particles on a swapping trajectory crosses the horizontal plane Az = in the 
direction corresponding to the counter-rotation of the vector connecting the particle cen- 
ters (cf. section I3.2[) . Since the counter-rotation ceases at the points where AU Z changes 
sign, the region of swapping trajectories is delimited by a set of critical trajectories for 
which AU Z — > when the particle pair approaches the plane Az = 0. (We note that for 
Az = the condition AU Z = implies that AU = because the horizontal components 
of the relative velocity (|3 .If) vanish because the symmetry of the problem.) 

For a two-wall system a typical graph showing the sign of the vertical relative velocity 
AU Z and the location of points where AU Z — in the now-vorticity plane Az = is 
presented in figure [H The set of points where AU Z vanishes consists of two subsets. One 
is the line Ax — in the region with no particle overlap, and the other is the circle 
Ap = A cl it (where Ap = Axe x + Aye y denotes the horizontal relative position vector, 
Ap = \Ap\, and A cr i t is the circle radius). The critical interparticle distance A cr j t does 
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Figure 11. Vertical extent Az of the upstream domains of swapping trajectories for particles 
moving symmetrically with respect to the midplane z — H/2 (as in the top panel of figure [10)) . 
plotted versus the channel width H. Solid line represents numerical results and dashed line the 
asymptotic behavior (|3.7p . The results in the inset are rescaled accordingly. 



not depend on the orientation of the particle pair in the plane Az — due to symmetry 
(cf. explanation in Appendix [BJ . 

By analyzing the direction of sphere motion, one can show that the area inside the 
circle of zeros Ap = A cr ;t corresponds to closed trajectories passing through the plane 
Az = 0. Therefore, only the part of the line Ax = outside this circle is associated 
with the critical trajectories that deli mit the swapping region. W e note that closed orbits 



egii 

also exist in unbounded shear flow (jBatchelor &: Greenl Il972t h In this case the whole 



surface Az = corresponds to particle pairs on closed trajectories (as in the top panel 
of figure [3]). 

For weakly confined systems the location of zeros of AU Z results from the balance 
between the ~ (Ap/d)~ A potential- flow contribution in equation (|3.3p and the wall re- 
flection of the stresslet. As already indicated in section [3~2l (cf. also Appendix [A)l . the 
leading-order behavior of the vertical component of the reflected flow in the region near 
the particle is 

Thus, for a fixed position z/H of the particle pair in the channel we find that 

A crit /d - (H/d) 3 / 5 , H/d > 1. (3.6) 

The dependence of the radius A cr jt on the channel width for a particle pair in the 
midplane z/H = | is plotted along with the asymptotic behavior (|3.6J) in the top panel 
of figure [9] The two bottom panels of figure [9] show the dependence of A cr i t on the 
position of the particle pair with respect to the channel walls for different H . The results 
represented in the left panel correspond to larger values of H/d, and they are rescaled 
to emphasize the asymptotic behavior (|3.6p . In the right panel the results are shown 
unsealed for moderate values of H/d. Note that for sufficiently small values of z/H the 
unsealed curves corresponding to different channel widths coincide, which indicates that 
the particle mobility is dominated by the single-wall contribution except for particles 
close to the center of the channel. 

An alternative way of characterizing the domain of swapping trajectories is to show the 
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Figure 12. Left panel represents domains in the upstream transverse coordinates (Ay, Az), 
corresponding to several values of the minimal gap e (as labeled) for a two-wall system with 
wall separation H/d = 5. Right panel shows the ratio between the frequency \ € of near contact 
particle encounters and the frequency Xawap of swapping-trajectory encounters in a dilute sus- 
pension under shear, versus the cutoff value of the minimal gap e for several wall separations 
H/d. All results are for trajectories symmetric with respect to the midplane z = H/2. 



corresponding upstream region of the transverse relative coordinates (Ay, Az), i.e., the 
region through which the swapping trajectories pass in the limit of infinite streamwisc 
particle separations Ax — > — oo. For trajectories that cross a given horizontal plane 
Az = during the particle motion, the border of the upstream region can be determined 
by integrating the trajectories backwards, starting near the loci of zeros of AU Z in the 
plane Az = 0. 

In figure [TU] the upstream boundaries of the swapping-trajectory regions are shown for 
several values of the wall separation H/ d. The top panel depicts these boundaries for the 
trajectories that are symmetric with respect to the midplane z = H/2 (strictly speaking, 
symmetric with respect to the axis defined by the intersection of the planes y = (y!+y 2 )/2 
and z = H/2). The bottom panels illustrate the dependence of the vertical extent of the 
upstream regions Az on the initial position z\ of the lower particle for particle pairs 
moving in the flow-gradient plane Ay = 0. In all the panels the swapping-trajectory 
domains correspond to the areas below the curves. 

By the scaling argument outlined in Appendix [C] one can show that the vertical extent 
of the upstream swapping-trajectory region scales as 

Az/d~ (H/d)- 1/2 , H/d>l, (3.7) 

which is confirmed by the results shown in figure 1111 The vertical axes of the plots 
shown in the top and left bottom panels of figure \W\ are scaled accordingly, to collapse 
the results for large values of H/d onto the asymptotic master curves. The results in 
the right bottom panel are shown unsealed to emphasize the near-wall behavior of the 
swapping trajectory domain. 

An interesting feature of the swapping-trajectory regions shown in the top panel of 
figure[TU]is that these domains are delimited by straight horizontal lines for Ay < Ay cr j t « 
A crit . The straight-line sections correspond to the trajectories that pass through the circle 
of zeros Ap = A cr ;t of the vertical relative velocity AU Z when the particle pair crosses 
the horizontal plane Az = 0. The remaining portions of the upstream borders of the 
swapping-trajectory regions correspond to the zeros of AU Z located at the axis Ax = 
in the plane Az = (cf. the geometry of zeros of AU Z depicted in figured]). The upstream 
coordinate Az in the region Ay < Ay cr jt is independent of the corresponding position 
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along the circle of zeros of AU Z for the same reason why A cr it does not depend on the 
orientation angle <j> (cf. Appendix |B|) . 

The results in the top panel of figure [TU1 indicate that the lateral extent Ay/H of 
the upstream swapping-trajectory region (normalized by the wall separation H) remains 
approximately constant when the wall separation H increases. However, the horizontal 
coordinate Ay CI - lt /H of the end point of the straight-line section decreases, as required 
by the scaling (|3.6[) . 

We conclude this section by discussing the distance of minimal interparticle approach 
for pairs of spheres on swapping and non-swapping trajectories. As we have already indi- 
cated, the interparticle gap (|3.2[) on trajectories of the swapping kind typically remains 
relatively large. For spheres that pass through a given horizontal plane Az — 0, the 
minimal interparticle separation on swapping trajectories coincides with the radius A cr it 
of the circle of zeros of the relative transverse velocity AU Z . According to the results 
shown in the top panel of figure [51 A cr i t is typically well above the particle diameter 
d (except for very small particle-wall separations). Thus, the minimal interparticle gap 
Emm = A cr ;t /d — 1 is usually O(l). In contrast, the minimal interparticle gap on non- 
swapping open trajectories may be smaller by several orders of magnitude. The sharp 
transition from the large gaps on the swapping trajectories to much smaller values of e 
on the nearby non-swapping ones is illustrated in the bottom panel of figure [3] 

The domains of non-swapping trajectories corresponding to specific values of the min- 
imal gap are depicted in the left panel of figure [12] for a system with a wall separation 
H/d = 5 and particles moving symmetrically with respect to the midplane z = H/2. 
The results indicate that the smallest values of the minimal gap are attained for near- 
critical non-swapping trajectories in the flow-gradient plane Ay — 0. According to our 
additional calculations the minimal gap decreases for stronger confinements, although it 
typically remains comparable to the distance of minimal appro ach emm = 4.2 x 10~ 5 for 
two spheres on open trajectories in the unbounded shear flow ([Arp fc Masonlll977l ). 

In a system of spheres interacting via a short-range repulsive potential, the near-contact 
particle encounters may result in cross-streamline particle displacements. However, for 
moderate-width channels, the upstream region corresponding to very small gaps is sig- 
nificantly smaller than the upstream region of swapping trajectories. 

Quantitatively, the rate \ °f particle encounters of a given type can be evaluated by 
integrating the upstream relative particle velocity AVoo over the appropriate upstream- 
coordinate region. The ratio Xe/xswap between the rates of near-contact and swapping 
collisions for trajectories symmetric with respect to the midplane z = H/2 is plotted in 
the right panel of figure [T^] versus the cutoff value of the minimal gap for channels of 
different width. The results indicate that x e /xswap <C 1 for sufficiently small values of 
the gap, especially at moderate confinements. 



4. Particle migration due to position-swapping binary encounters 

As shown above, the position-swapping particle encounters result in cross-streamline 
particle migration in suspensions confined by planar (or nearly planar) walls. Although 
this migration mechanism may be dominant in dilute suspensions under moderate-confinement 
conditions, it has never been described. In this section we present our quantitative pre- 
dictions for migration of spherical particles in parallel wall channels. 

4.1. Population-balance simulations 
In order to illustrate the role of particle swapping in suspension flows, we have performed 
a population-balance simulation of a mixture of two species of equal-size spheres in a 
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Figure 13. Particle migration in a dilute binary suspension of equal-size spheres undergoing 
shear flow in a parallel- wall channel with wall separation H/d — 5. Top panels show the initial 
state of the system, and bottom panels represent the state after 128 swapping binary encounters 
per particle. The left panels depict vertical positions z/d for randomly selected 200 particles of 
each species (versus the particle index k), and right panels show the normalized particle density 
n across the channel. The curves in right-bottom panel represent the 3 rd order polynomial fit 
to the histogram of the particle distribution (taking the system symmetry into account) 
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Figure 14. Swapping-trajectory contribution (|4.6[) to the hydrodynamic self-diffusivity coeffi- 
cient. Left panel shows Dg Wap versus distance z/d from the lower wall for different values of the 
channel width H/d (as labeled). Right panel shows _Ds Wap at the center of the channel versus 
H/d. Dashed line represents the asymptotic behavior (|4.7jl . and the inset shows the results in 
the corresponding rescaled variables. 
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parallel-wall channel of width H/d = 5. Initially, the particles of species one are located 
in the upper portion of the channel, z > H/2, and the particles of species two in the lower 
portion z < H/2. In each part of the channel the particles are distributed randomly. Since 
the system is translationally invariant in the transverse directions x and y, the particle 
distribution is characterized only by the vertical coordinate z. 

To mimic the evolution of a dilute suspension, the particle population is updated 
using a Boltzmann-Monte Carlo technique. Accordingly, a pair of particles is chosen 
randomly from the ensemble representing the system, with a probability proportional 
to the upstream relative velocity At/oo- Then, a two-particle evolution run is performed 
with a large initial offset in the streamwise direction |Ax| ^> d and the transverse offset 
Ay chosen randomly from the uniform probability distribution truncated at Ay max w 
8d. If the pair evolution results in position swapping, the particle ensemble is updated 
accordingly. 

Sample results of our population-balance simulations are presented in figure 1131 The 
calculations were performed using a set of 600 particles, and the average over seven 
simulation runs was taken to reduce statistical fluctuations. The top panels show the 
initial random particle distribution, and the bottom panels the distribution after 128 
swapping collisions per particle. In the left panels the distribution is depicted by plotting 
the vertical positions of 400 randomly chosen particles versus particle index. The right 
panels represent the corresponding particle density (normalized to unity) as a function 
of the position z across the channel. The results indicate that the two particle species 
slowly mix together due to the particle-swapping phenomenon. Without particle swap- 
ping, all spheres would preserve their initial vertical positions (unless non-hydrodynamic 
interparticle forces are present). 



4.2. Self- diffusion coefficient 

Assuming that the wall separation H is much larger than the swapping range Az, a 
sequence of uncorrelated particle displacements due to binary encounters in suspension 
flow can be approximately described as a diffusion process. In this section we focus on the 
transverse diffusivity of a tagged particle in a dilute suspension of mechanically identical 
spheres. 

As shown by Ida Cunha fc Hindi (1996) in their paper on the shear-induced sel f- 



diffusivity in a dilute suspension of rough spheres (also see IZarraga k Leightonll200lh . 



the transverse (cross-streamline) self-diffusion coefficient D s , can be evaluated from the 
relation 



\n J [(AZfAU^dAydAz, (4.1) 

where the integration is over the upstream region of transverse relative coordinates cor- 
responding to particle encounters of a given type (e.g., swapping or direct particle colli- 
sions), n is the particle number density and AZ is the cross-streamline displacement of 
the particle during a binary-encounter event. 

Since initially the particles arc well separated in the streamline direction, the relative 
upstream velocity AU<x, = Uoo{z2) — Uoo{zi) is the difference of the velocities Uoc(z) of 
individual, non-interacting spheres at the positions z\ and Z2 with respect to the channel 
walls. For particles far from the walls we simply have AUoa = 7A2:. More generally, 
assuming that the initial particle offset Az is sufficiently small, the relative particle 
velocity AUoo can be expressed as 

AUoo = jaAz, (4.2) 
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where 

*{z)=T 1 ^ (4-3) 
dz 

is the dimensionless derivative of the velocity of an individual sphere with respect to it s 



transverse position in the channel. According to our results ()Zurita-Gotor et aZ.1 120071 ) . 
the velocity of a single sphere undergoing shear flow in a parallel- wall channel differs very 
little from the local fluid velocity ()2.2p , except when the particle is in the immediate wall 
proximity. Therefore, the approximation a = 1 is often sufficient. 

Particles on swapping trajectories exchange their vertical positions, which implies that 
the particle displacement AZ equals the upstream offset Az of the particles in the in- 
teracting pair. Using the above result and relation (14. 2|) . the integral with respect to 
Az in equation (|4.ip can be explicitly performed over the swapping trajectory region. 
The effective self-diffusivity can thus be determined from the resulting one-dimensional 
integral 

/oo 
A^ ax dAy, (4.4) 
-oo 

where Az max = Az max (Ay) is the upstream vertical particle offset on the trajectory 
that delimits the region of swapping trajectories (cf. the contour plots in the top panel 
of figure [TO]) . Rephrased in terms of the volume fraction <f> = ^Trd 3 n and dimensionless 
particle offsets Ay = Ay/d and Az = Az/d, equation (|4.4|) yields 

D s s wap = ±jcj)d 2 D s s wap , (4.5) 

where 

/•OO 

5 r a P = 67r -i / A3£ ax dAy. (4.6) 



The dimensionless self-diffusion coefficient £)^ wa P depends on the position z of the 
particle pair across the channel. Since the problem is non-local, we assign z to be the 
position of the particle pair when it crosses the horizontal plane Az = 0. 

Our numerical results for the swapping contribution (14. 6|) to the effective self-diffusion 
coefficient are presented in figure Q31 The left panel shows D* wa -P versus z for several 
values of channel widths, and the right panel shows the self-diffusivity at the channel 
center as a function of H . The results of our calculations indicate that for H/d < 15 the 
self-diffusion coefficient D;! wap achieves its maximal value at the channel center, where 
the swapping effect is the strongest due to the superposition of the flow reflected from 
two walls. In contrast, for H/d > 15 the maximum occurs at the distance z/d « 2 from 
each wall. The results also show that particle migration is suppressed near the walls, 
which is consistent with the observation that the vertical extent of the particle-swapping 
domain is small near the walls (cf . the bottom panels of figure [TO]) . The effect of the 
decreased self-diffusivity in the near-wall regions is reflected in the shapes of the density 
profiles depicted in figure fT3l 

The behavior of £)g Wap for H/d ^> 1 can be determined by combining relations f|3. T|) 
and (|4.6j) and recalling that the lateral extent Ay of the swapping region scales with H 
for the large channel width (as depicted in the top panel of figure [T0|) . Accordingly, we 
obtain the relation 

= qj^/H, (4.7) 

where q depends on the position of the particle pair across the channel z/H . For particles 
in the midplane z/H — ^ we find q ~ 0.55, according to the results shown in the inset 
of the right panel of figure HU 
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4.3. Comparison with experiment by Zarraga and Leighton 

While the mechanism for particle migration due to the swapping traj ectories has not been 
propo sed so far, the effect of the swapping motions was observed when lZarraga fc Leighton 
( 20021 ) reported unusually large values of the hydrodynamic self-diffusivity for a dilute 
suspension of spheres undergoing shear flow in a Couette device. Their measurements 
gave the result D s = 3.6 x 10~ 2 70a 2 (where a = is the particle radius) for the low- 
density self-diffusion coefficient, whereas the estimate of the contribution due to particle 
roughness was at least four times smaller: it ranged from _D™ ugh = 3.0 x 10~ 3 70a 2 to 
8.4 x lO~ 3 70a 2 , depending on the ass umed roughness amplitude. This discrepancy has 
been unaccounted for, and the result of lZarraga fc Leightonl (|2002T ) experiment has been 
puzzling. 

By analyzing the effect of the channel walls on particle motion we provide a simple ex- 
planat ion for the anomalous value of the self-diffusion coefficient. The lZarraga fc Leighton 
(2002) measurements were performed in a channel with wall separation H/d = 20. For 
this geometry we find that the swapping contribution to the self-diffusion coefficient in the 
center of a channel is £)| wa P — 2.4 x 10~ 2 70a 2 . Moreover, in the whole central region the 
diffusivity £)| wa P only weakly depends on the position z, according to the results shown in 
figurc [T4l Assuming t he up per limit of the roughness parameter e r = 1.8 x 10~ 2 quoted by 
IZarraga fe Leightonl (|2002f) . we find L>™ugh = 6.7 x 10" 3 70a 2 (note that £>™ u s h is slightly 
smaller than the corresponding value in infinite space because the swapping mechanism 
reduces the rate of direct particle collisions). Combining the swapping and roughness 
contributions, we obtain D s = 3.1 x lO~ 3 70a 2 for the total self-diffusion coefficient. Our 
result agrees with the experimental value with accuracy of 15 %. 

In this paper we focus on systems of equal-size spheres but wall- induced cross-streamline 
particle migration also occurs in bidisperse or polydisperse suspensions. We note that 
in such suspensions, migrati on resulting from particle rough ness is significantly smaller 
than in monodisperse ones (|Zarraga fc Leightonl [200ll l2002h . Thus in polydisperse sys- 
tems particle swapping constitutes the dominant particle migration mechanism even if 
there is significant particle roughness. 

Moreover, cross-streamline displacements of two different-size particles undergoing a 
binary swapping encounter differ because there is no fore-aft symmetry. It follows that 
the total suspension density is affected by the swapping trajectories, in addition to the 
migration of individual particle species. In contrast, in monodisperse suspensions the 
total density cannot be altered by position-swapping binary encounters. 



5. Conclusions 

The key finding of our study is that confining walls can qualitatively change the topol- 
ogy of binary encounters of spherical particles in suspension flows. Specifically, we have 
identified a new class of binary trajectories that result in cross-streamline particle mi- 
gration in a wall-bounded shear flow. Equal-size spheres on such trajectories do not pass 
each other (as in an unbounded system) but, instead, they exchange their transverse 
positions. While our explicit calculations arc limited to equal-size spherical particles in 
shear flow bounded by parallel planar walls, we expect that a similar effect also exists 
for pressure-driven flows, particles of different sizes, and other wall geometries. 

We have shown that the cross-streamline particle motion is driven by the wall reflection 
of the scattered flow produced by the particles. Namely, the reflection of the flow produced 
by one of the spheres pushes the other approaching sphere across the streamlines of the 
external flow towards the fluid region moving in the opposite direction. Due to this 
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mutual interaction, the particles turn around and return to infinity without passing each 
other - they exchange their vertical positions instead. 

The significance of our results stems from the fact that the wall-induced transverse 
particle displacements associated with position-swapping binary encounters constitute 
the sole mechanism for cross-streamline migration of spherical non-Brownian particles 
at low suspension concentrati ons (if there are no n o n-hyd rodynamic forces) . Such migra- 
tion was indeed observed by IZarraga &: Leighton ( 20021 ) who found anomalously large 
self-diffusion coefficient in a confined suspension of spheres. Up till now, however, their 
measurements have remained unexplained. In this paper we show that the swapping 
mechanism quantitatively explains the unexpected results of their experiments. 

A sequence of uncorrelated binary position-swapping particle encounters in a monodis- 
perse suspension causes a macroscopic self-diffusion process. In multi-component systems 
a sequence of such encounters results in mutual diffusivity. Moreover, the swapping tra- 
jectories produce migration not only in the velocity-gradient direction (considered herein) 
but also in the direction of fluid vorticity. 

In addition to our two-particle results, we have found that in the wall presence there 
is a domain of recirculating streamlines around a single sphere in shear flow. Such recir- 
culating streamlines result in fluid mixing in suspension flows, so this phenomenon may 
have important microfluidic applications. Furthermore, the trajectory reversal associated 
with particle swapping behavior may prevent near-contact particle encounters by causing 
the particles to separate before arriving into a near-contact configuration. Such a hydro- 
dynamic shielding effect may be especially important for a chain of particles with slightly 
different cross-streamline positions in a microfluidic channel. Thus, proper understanding 
of the position-swapping mechanism may contribute towards finding better methods of 
controlling particle motions in microfluidic devices. 

Other consequences of the particle-swapping mechanism will be described in our forth- 
coming publications. In particular, we will show that binary particle collisions in a dilute 

suspension under shear may produce a layered particle distribution. Our preliminary re- 

suits concerning such a layering process have already been presented ( Zurita-Gotor. Blawzdziewicz fc Wainrvbl 
2005h . 
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Appendix A. Flow reflected from rigid wall 

According to the lLorenta (jl907T ) expression, the wall reflection 6y* of the flow field (13.31) 

in a r igid wall at z = can be expressed by the following formulas ( Bhattacharva fc Blawzdziewicz! 
20021) . 



K = P(Ro + ziRi + zfR 2 ) • <5 Vl 



where 



Rn 



2zVe z + z 2 V 2 I, 



(Al) 
(A 2) 
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Ri = -2Ve z + 2zV 2 I, (A 3) 

R 2 = V 2 I. (A 4) 

In the above relations, z — z — z\ is the coordinate relative to the source position, I is 
the identity tensor, I z = I — 2e z e z , and P is the reflection operator with respect to the 
plane z = 0, 

[Pw](x,y,z) = I z ■ w(x,y, -z). (A 5) 

By inserting f|3. 3|) into (|A 1[) and evaluating the order of magnitude of different terms, 
we find that the overall decay of the reflected flow Sv* with the wall-particle distance z\ 
is 0((d/ zi) 2 ). Equations (|3.3|) and (|A ip also indicate that the vertical component 8vl z 
is an odd function of the variable x. Therefore, for x/ ' z\ < 1 we find 




(A 6) 



For a given polar angle <f>, equation (|3 . 5[) is obtained from (|A 6[) and an assumption 
that z\ ~ H. To obtain relation (|3.6|) . we use (|3.3| . (|3.5p . and the observation that 
the dependence of both 5v\ z and 5v* z on the polar angle <f> is cos (which follows from 
symmetry, as explained in Appendix [Bj . 



Appendix B. Location of zeros of AU Z in the plane Az = 

The critical interparticle distance A cr ;t does not depend on the orientation of the 
particle pair in the plane Az = 0, which can be shown by decomposing the external flow 
(|2.2|) into the longitudinal component (along the line connecting particle centers) and 
transverse component (normal to this line), 

y ext =v ext +v ext_ (B 1) 

By symmetry, the transverse velocity v^* does not produce any vertical particle motion. 
Therefore, only the longitudinal problem defines the location of zeroes of AU Z . Moreover, 
the location of points AU Z = does not vary with the magnitude of the projection vjj xt . 
Thus the distance A cr ;t is independent of the orientation of the vector Ap. 

By a similar argument one can demonstrate that the upstream swapping-trajectory 
domain in the variables (Ay, Az) is delimited by a straight line section Az = const for 
Ay < Ay cr jt (as shown in the top panel of figure [T0|) . The straight-line section corre- 
sponds to trajectories that pass through the circle of zeros of AU Z . Since the transverse 
component Vjf* of the external flow does not produce any vertical particle motion, for 
all initial orientations of the particle pair on the circle Ap = A cr ; t the vertical evolution 
of the particles is thus the same, except for a varying time scale. 

Using the flow decomposition (|B ip one can also show that for a given Ap, the vertical 
velocity AU Z varies as cos <j> with the polar angle 4> of the relative-position vector Ap 
(AU Z varies with <j> in the same way as the flow component vjj xt ). The plots of AU Z vs 
Ax shown in figured] are thus sufficient to determine AU Z for all orientations of a particle 
pair in a given plane z — const. 

We would also like to mention that our numerical results indicate that additional 
changes of sign in the relative vertical velocity AU Z may occur at large interparticle 
distances. This would lead to more complex particle recirculation patterns than those 
described herein. However, the magnitude of AU Z a t large distances is extreme ly small 
due to the exponential decay of the vertical velocity ( Bhattacharva et a?Jl2006ah . There- 



fore, domains of trajectories with different topology are limited to very small regions in 
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the configurational space. We note tha t multiple c hange s in sign of the velocity field in a 
parallel-wall channel were described bv lHackbornl (|l990j ) for a flow produced by a vertical 
rotlet. 



Appendix C. Proof of relation (13. 7|) 

As we have shown in our previous publications ( Bhattacharva et al. 200 6a. b). the 
vertical velocity components for two spheres in Stokes flow in a parallel-wall channel 
decay exponentially on the lengthscale H . Thus, the vertical offset of the particles along 
the pair trajectory evolves only for Ax < H . The relative lateral velocity of the particles 
in shear flow (|2.2[) is AU X ~ jAz, where Az is a typical vertical offset along the trajectory. 
Therefore, 

to-j^H/Az (CI) 

is the timescale for the evolution of the vertical particle positions. During the time interval 
(jC II) the vertical offset Az may attain the magnitude 

Az~t AU Zl (C2) 

where AU Z is the typical value of the vertical component of the relative particle velocity. 

According to equation f|3. 3|) . the perturbation velocity field Sv\ produced by a spher- 
ical particle of diameter d in an unbounded shear flow with the shear rate 7 scales as 
Svi ~ jd(d/f) 2 , where r is the distance from the particle. In a wall-bounded system, 
this perturbation flow interacts with the walls at the distance f ~ H and then acts on 
the second particle at a distance Ax < H, producing its vertical motion. Therefore, the 
relative vertical velocity of the particles has the magnitude 

AU z ~jd{d/H) 2 . (C3) 

Relation (|3.7|) is obtained by combining equations (jC 1|) - (|C 3jl . 
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